#!/bin/bash
#
# Script to calculate the surface energy for Cu(001)
#
# gamma = (E_slab_unrelaxed - N * E_bulk_per_atom) / (2*A) + (E_slab_relaxed - E_slab_unrelaxed) / A
#

E0=-3.73075857218744  # Energy for bulk FCC Cu, eV/atom
#
# OUTCAR.01 and OUTCAR.03 are needed.
if [ ! -f "OUTCAR.01" ]; then
   echo "ERROR: File OUTCAR.01 not found, calculation incomplete!"
   exit
fi
if [ ! -f "OUTCAR.03" ]; then
   echo "ERROR: File OUTCAR.03 not found, calculation incomplete!"
   exit
fi
# get surface area info
a0=`head -2 POSCAR|tail -1|awk '{print $1}'`
lx=`head -3 POSCAR|tail -1|awk -v s=${a0} '{print $1*s}'`
ly=`head -4 POSCAR|tail -1|awk -v s=${a0} '{print $2*s}'`
area=`echo ${lx} ${ly}|awk '{print $1*$2}'`
#
# get number of atoms info
natom=`grep 'NIONS' OUTCAR.01|tail -1|awk '{print $12}'`

# get energy difference between ideal surface and bulk material
Eunrel=`grep 'free  energy' OUTCAR.01|tail -1|awk '{print $5}'`
DelE=`echo ${Eunrel} ${E0} ${natom}|awk '{print $1 - $2*$3}'`

# get energy difference between the relaxed and ideal surface
Erelax=`grep 'free  energy' OUTCAR.03|tail -1|awk '{print $5}'`
DErel=`echo ${Erelax} ${Eunrel}|awk '{print $1 - $2}'`

# surface energy
gunrel=`echo ${DelE} ${area}|awk '{print $1/(2*$2)}'`
gamma=`echo ${DelE} ${DErel} ${area}|awk '{print $1/(2*$3) + $2/$3}'`
echo ""
echo "The energy of the as-cut  slab is: ${Eunrel} eV."
echo "The energy of the relaxed slab is: ${Erelax} eV."
echo "The surface area of the   slab is: ${area} Ang^2."
echo "The unrelaxed  surface  energy is: ${gunrel} eV/Ang^2."
echo "The relaxed    surface  energy is: ${gamma} eV/Ang^2."
echo ""
exit 0
